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The notion of a dynamical quantum phase transition (DQPT) was recently introduced in [Hcyl 
et al., Phys. Rev. Lett. 110 , 135704 (2013)] as the non-analytic behavior of the Loschmidt echo 
at critical times in the thermodynamic limit. In this work the quench dynamics in the ground 
state sector of the two-dimensional Kitaev honeycomb model are studied regarding the occurrence 
of DQPTs. For general two-dimensional systems of BCS-type it is demonstrated how the zeros 
of the Loschmidt echo coalesce to areas in the thermodynamic limit, implying that DQPTs occur 
as discontinuities in the second derivative. In the Kitaev honeycomb model DQPTs appear after 
quenches across a phase boundary or within the massless phase. In the Id limit of the Kitaev 
honeycomb model it becomes clear that the discontinuity in the higher derivative is intimately 
related to the higher dimensionality of the non-degenerate model. Moreover, there is a strong 
connection between the stationary value of the rate function of the Loschmidt echo after long times 
and the occurrence of DQPTs in this model. 


I. INTRODUCTION 

Recent advances in experimental techniques allow to 
realize closed quantum systems with cold atomic gases 
in optical traps A— These setups are precisely controllable 
and the unitary time evolution of the systems can be re¬ 
solved such that the dynamics are experimentally acces¬ 
sible under well-known conditions. Motivated by the new 
experimental possibilities a lot of theoretical research on 
the non-equilibrium dynamics of quantum systems has 
been conducted in the past years. In these theoretical 
investigations a common protocol for driving a system 
out of equilibrium is called quantum quench. Considering 
a parametrised Hamiltonian H(a), where the parameter 
typically corresponds to some external field strength in 
the experimental setup, the system is initially assumed 
to be in equilibrium with regard to some value a, of the 
parameter. Then, the parameter is suddenly quenched 
to a different final value a/ driving the system out of 
equilibrium and inducing a non-trivial time evolution. 

Studying the quench dynamics of a quantum many- 
body system, Heyl et al4 pointed out the close formal 
similarity of the canonical partition function of an equi¬ 
librium system, Z(f3) = tr and the return ampli¬ 

tude 

Q{t) = (iPi\e- lHt Wi) (1) 

of a time-evolved state, suggesting the possibility of crit¬ 
ical behavior in the time evolution in analogy to equilib¬ 
rium phase transitions. It is known that in the thermo¬ 
dynamic limit the zeros of a partition function coalesce 
to lines in the complex temperature plane and the equi¬ 
librium phase transition is marked by the intersection of 
the zero-line with the real temperature axis^ Heyl et al. 
found that in the case of the transverse field Ising model 
the boundary partition function 

Z{z) = (i/)i\e~ zH ( 2 ) 

has zeros in the complex time plane, which accordingly 
coalesce to lines in the thermodynamic limit. These lines 


cross the real time axis after quenching the external field 
across the quantum critical point inducing nonanalytici- 
ties in the rate function of the Loschmidt echo 

r(t) = - lim 2 \n\(ipi\e- im \ipi )\ 2 

iV—^OO iV 

= -Jim K7 ln£ W ( 3 ) 

at equidistant critical times f*. Heyl et al. denote this 
non-analytic behavior at critical times in the thermody¬ 
namic limit as a dynamical quantum phase transition 
(DQPT). They showed that in experiment the DQPT 
would be observable by measuring the work distribution 
function of a double quench; in particular, the Loschmidt 
echo C{t) = \(4>i\e~ lHt \ifi )\ 2 is the probability of perform¬ 
ing no work. 

These findings triggered further work aiming at a bet¬ 
ter understanding of the phenomenon. By considering an 
additional integrability-breaking interaction in the trans¬ 
verse held Ising chain it was demonstrated that DQPTs 
are not a peculiarity specific to integrable models, but 
are stable against some non-integrable perturbations;^ 
Moreover, the signature of DQPTs was found in higher¬ 
dimensional systems, namely, in two-dimensional topo¬ 
logical insulators^ and effectively infinite dimensions us¬ 
ing DMFT£ It was observed in various cases that DQPTs 
are not necessarily connected to quenching across a quan¬ 
tum critical pointer— however, there seems to be a 
strong connection to topological phase transitions 
Canovi et al4 detected coexisting solutions for so called 
generalized expectation values in post quench dynamics 
and, therefore, they introduced the notion of a first order 
dynamical phase transition. This could be a way to clas¬ 
sify dynamical phase transitions. Furthermore, a close 
connection between the analytic behavior of r(t ) in the 
complex plane and its long time limit is conjectured; 13 

In this work we study quench dynamics in the Ki¬ 
taev honeycomb model^ regarding dynamical quantum 
phase transitions. The model features a rich phase di¬ 
agram comprising an extended gapless phase, anyonic 
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excitations, and topological order. Moreover, it is a 
rare example of a Jordan-Wigner-solvable model in two 
dimensions . 15 ^ 16 As such it has been studied extensively 
under various aspects. In this paper we restrict the dis¬ 
cussion to the dynamics in the ground state sector. 

The dynamics of gapped two-dimensional two-band 
systems were already studied by Vajna and Dora with 
focus on a connection between DQPTs and topological 
phases. 7 In the presence of a magnetic field the Kitaev 
model acquires topological order and becomes a system 
of the same family, albeit being a spin model. In that 
case we find the behavior in accordance with their re¬ 
sults, namely, DQPTs occur after quenches across the 
boundary between phases with different Chern number. 
However, the focus of this work lies on quenching be¬ 
tween the topologically trivial phases in the absence of a 
magnetic field. Similar to their results we find DQPTs 
as discontinuities in the second derivative, which is in¬ 
herent to the higher dimensionality of the system, and 
we elaborate on the relevance of the complex zeros of the 
dynamical partition function in this context. Moreover, 
we discuss a remarkable observation regarding the long 
time behavior of the Loschmidt echo. If no DQPTs oc¬ 
cur in the post-quencli dynamics, then, although the ap¬ 
proached stationary state is always an excited state, the 
long time limit of the Loschmidt echo is given by the fi¬ 
delity, i.e., the overlap of the initial state with the ground 
state of the quenched Hamiltonian. This, however, does 
not hold if the dynamics exhibit DQPTs. 

The rest of this paper is organized as follows: in sec- 
tion[TT]the way of solving the model using Jordan-Wigner 
transformation is sketched and the phase diagram is in¬ 
troduced. Furthermore, the expressions for the dynami¬ 
cal free energy is derived. In section Hill the zeros of the 
dynamical partition function in the complex time plane 
are treated assuming a general BCS-type Hamiltonian, 
yielding the criteria for the occurrence of DQPTs and 
the order of the corresponding nonanalyticity. Finally, 
the zeros of the partition function and the real time evo¬ 
lution are studied explicitly for the Kitaev model in sec¬ 
tion IIV1 and two interesting limits are taken into account 
as well as ramping as an alternative protocol and the 
quenching with an additional magnetic field. 


II. THE KITAEV HONEYCOMB MODEL 
A. The model 

The Kitaev honeycomb model is defined by the Hamil¬ 
tonian 


H(f) = - E E J " a o :<J k , ( 4 ) 

aE:{x,y,z} a-links 


which describes a spin-1/2 system with the spins located 
on the vertices (labeled by j, k) of a honeycomb lattice as 
depicted in Fig. [TJ 1 - In this paper we assume the lattice 



FIG. 1. Lattice of the Kitaev honeycomb model given by eq. 
©• Spin-1/2 degrees of freedom are sitting on the vertices of a 
honeycomb lattice. The nearest neighbor interaction depends 
on the link type ( x , y, or z). 


spacing to equal unity. It has been showed that for the 
above Hamiltonian one can find a Jordan-Wigner con¬ 
tour, which after identifying a conserved Z 2 operatoi— 
and switching to momentum space yields a BCS-type 
Hamiltonian 






(^-cLjdis) 




( 5 ) 


with 

e^r(t7) — 2(^I Z T J x cos(/c^) T Jy cos(/cy)) , 

Ag( J) = 2 (J x sin (k x ) + J y sin(£: y )) . (6) 


This Hamiltonian can be diagonalised by a Bogoliubov 
transformation 



( u k^) v k(J) \ f ( ¥ \ 

\~ v k( J )* \ d -k) ’ 


( 7 ) 


where 


u k(J) 


i( i+ m 

2 V e-aj) 


v k(J) = sg n( A fc(y)) 


1 ( x yj) 

2 l EAJ) 


( 8 ) 


(see appendixlAlfor details). Plugging the transformation 
into eq. (jn]) yields the diagonal Hamiltonian 

HW = E ( a h - «-tol E ) (9) 

keK 


with spectrum 

E^J) = ^/e^jy + A^jy . ( 10 ) 
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FIG. 2. Phase diagram of the Kitaev model in the J x + J y + 
J z = 1 plane, where J a > 0. The gapless phase B is sur¬ 
rounded by three gapped phases A a , where J a > Jp + J 1 . 


The Hamiltonian splits into a sum over fc-sectors, i.e., a 
sum over k £ K , where K is a subset of the Brillouin- 
zone such that Vfc £ K , —k £ K. A possible choice is 
one half of the Brillouin-zone, e.g., all k with k x > 

The spectrum has roots at 

( Jl- J 2 - J 2 \ 
k x = ± arccos —- 1 , 

^ 2 JxJz ) ’ 

( J 2 -Jl- Jl \ 

k y = T arccos ( -- j (11) 

if |J a | < \Jp\ + |J 7 |, where is any permutation 

of (x,y,z). We will in the following only consider non¬ 
negative J a on the J x + J y + J z = 1 plane. In this sec¬ 
tion the above result on the gappedness of the spectrum 
corresponds to a phase diagram as depicted in Fig. [2] 
The gapless phase B at the center of the diagram is sur¬ 
rounded by three distinct gapped phases^ A x , A y and 
A z , where J x Jy T J -, Jy J x T J z or J z J y T J x , 
respectively. 

In the presence of a magnetic field h the spin Hamil¬ 
tonian 0) becomes 

H(J,h) = - ]T j £ JacrX + E^I 

<y.£{x,y,z} \ a-links j J 

(12) 

and the additional term opens a gap also in the H-phase. 
Moreover, the H-phase becomes topologically non-trivial 
with Chern-number v = ±1, whereas the A a phases re¬ 
main trivial with v = 0. 14 At J x = J y = J z = J there 
exists a diagonal form of the Hamiltonian even with non¬ 
zero magnetic field" 14 and the spectrum reads 

W h) = ^(J,«) 2 + |A S (J,«)| 2 (13) 

with ee( J, n) = eAJ), J = ( J , J, J) T and 
Afe( Ji K ) =A j:( J) 

+ 4 in (sin(fc x ) — sin(fc y ) + sin [k y — k x )) , 

(14) 


where n ~ x j \ hz . Through a Bogoliubov transforma¬ 
tion (see appendix [All this maps to 


H(J, K )=J2 


keK 


e ^ J ' ^ (didt -d - 

2 V k k ~ k -k 


A fc( J ’ K ) A Tf ' A k( J ’ K ^ d ^ 
2 —k k 


- did T 

2 k -k 


(15) 


This case will be studied in section IIVFI Before, we will 
stick to the case without magnetic field, i.e., real valued 


££ and Ag. 


B. Post-quench dynamics 


In order to study the dynamics in the Kitaev honey¬ 
comb model we will consider situations where the system 
is initially, at t < 0, prepared in the ground state of 
Jf(Jo), i.e., Ff(J 0 )|^>j} = Ecsl^i)- In terms of the free 
fermion degrees of freedom the initial state is the vac¬ 
uum: al°\ipi) = a^°|0; Jo) =0. At t = 0 the parameter 
is quenched to its final value J\, such that for t > 0 the 
time-evolved state is \ip(t)) = Making use 

of the Bogoliubov transformation the initial state can be 
expressed in terms of the final free fermions, which diag- 
onalise H(Ji): 


m = A/-” 1 n (i + %(Jo, fjafa&A |0; £) (16) 

fceic 


Here, 

E%{Jo, Ji) 


U^(JoJi) 

M fc( J o)ng(Ji) ~ ng(Ji)ng(J 0 ) 

u k(Ji) u k(Jo) + ^(Jl)^(Jo) 


(17) 


and the normalization constant 

j^= n i i + B k( j, o,Ji) 2 ) as) 

k£K 


were introduced. A more detailed derivation is given in 
appendix [A] 

For the sake of brevity and lucidity we will in the fol¬ 
lowing refrain from dragging along the dependencies on 
Jo and Ji explicitly, i.e., identify Bg = B^Jq, Ji) and 

a k= a t- 

Using m to compute the dynamical partition func¬ 
tion we obtain 


Z(z) = (ipi\e~ zH \ipi) 


= n 


1 + Ble~ 2E k( f ^ 

_ k _ 

l + Bl 


(19) 


( 20 ) 
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The dynamical partition function has large deviation 
form Z(z) ~ e~ N f( z \ where N is the system size. Thus, 
in the thermodynamic limit only the rate function, or 
dynamical free energy density, 

f( z ) = ~ I™ T j M z i z )) (21) 

Af->oo 7V 

is well defined. 


III. ZEROS OF THE PARTITION FUNCTION 
AND CRITICAL TIMES 

A. General aspects 



FIG. 3. Schematic picture of the surface separating two densi¬ 
ties of zeros in the complex plane and the relevant coordinate 
frames for determining the behavior of 4>{z) (eq. (1271) 1 along 
the real time axis. 


From the study of equilibrium phase transitions it is 
known that a very insightful approach is to consider the 
zeros of the partition function in the complex tempera¬ 
ture or complex magnetization plane, respectively . 4 ' 19 ! 20 
In the thermodynamic limit the zeros of the partition 
function coalesce to lines or areas in the complex plane, 
which mark the critical points when approaching the real 
temperature (magnetization) axis. Analogous reasoning 
has proven useful in the study of dynamical quantum 
phase transitions—— 

In particular, an interesting analogy to electrodynam¬ 
ics allows to characterize the dynamical phase transition 
through the density of zeros of the dynamical partition 
function in the complex time plane. The starting point 
is the observation that the dynamical partition function 
m is an entire function of z and can as such, according 
to the Weierstrass factorization theorem, be written as 

Z(z) = e h W Ilf 1 --) ’ ( 22 ) 

jeJ V ZjJ 

where J is some discrete index set, Zj £ C are the zeros, 
and h(z) is an entire function— With this, the dynamical 
free energy density reads 


/(*) 


lim 

N—too 


l 

N 


H z ) + ln ( 1 

j£J V 



(23) 


From this expression it becomes clear, that any non- 
analytic behavior of the dynamical free energy can only 
occur at or in the vicinity of the zeros of the dynamical 
partition function Zj. Since we are interested in nonana- 
lyticities, we will in the following ignore the contribution 
of h(z) and only consider the singular part 

f s (z) = - lim -T7 ln f 1 “ —) ( 24 ) 

iV—>-oo iV z —' V Zj J 

j&J v J/ 


In the thermodynamic limit the sum becomes an integral 
over some continuous variable x £ X, where X C R” is 
a region corresponding to the previously used index set 


J, and the zeros become a function of this variable z(x ), 
such that 


rw= -X dxh, ( 1 - ife) ■ (25) 

A transformation of the integration variable yields 

f s {z) — ~ f dzp(z) In (l - , (26) 

Jz(X) v z) 

where the Jacobian determinant p(z) can be interpreted 
as the density of zeros in the complex plane—I Moreover, 
setting p{z) = 0 for z z(X) allows to extend the inte¬ 
gration domain to the full complex plane. We will now 
discuss the real part 


<fi(z) = Re [ f s (z )] = [ d zp(z) ln 

Jc 


1 - Z - 

z 


(27) 


We will later see that the Loschmidt echo on the real 
time axis is directly given by <j>(t). For z = x + iy with 
x,y £ R ln \z\ is the Green’s function of the Laplacian 


A 2 d = 


dx' 2 


dy 


-, i-e., 


A -2d4>(z) = -27T p(z) . 


(28) 


In other words, the real part of the dynamical free energy 
density can be interpreted as the electrostatic potential 
cf>(z) produced by a charge density p(z) in two dimensions 
and the question of the behavior of the free energy at 
critical points becomes the question of the behavior of 
the electrostatic potential at surfaces. If the zeros form 
lines in the complex plane, this allows to deduce the order 
of the phase transition directly from the density of zeros 
at or in the vicinity of the physically relevant z> 

Although the zeros can form areas in the complex plane 
these areas do not cover the physical axis in the case of 
equilibrium phase transitions. However, as we will see in 
the following section, this is possible in the case of dy¬ 
namical phase transitions. Thus, consider the situation 
as depicted in Fig. [3] The density of zeros is given by 
Pi(z) and p 2 (z) in area I and II , respectively, and at the 
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boundary there is a discontinuous change in the density 
of zeros. Assume the electric potentials <j>i(z), i = 1,2, 
solve the Laplace equation (l28l) with the corresponding 
density Pi(z). With this a global solution is 




f<M~) ,z€l 
\ (z) ,z G II 


(29) 


and we demand continuity of <fi(z) at the boundary. Let 
us now focus on the behavior of (f>{z) at the intersection of 
the boundary with the real time axis. If we choose the x- 
y-coordinate frame as indicated in Fig. [51 the continuity 
of <j>(z) implies that on the boundary 

■^Mz) --^Mz) = 0 ■ (30) 


Since we are interested in the behavior of in the real 
time axis, we transform to '-coordinates, 


t = 


cos a 

y' = y ■ 

With the Laplace equation 

d 2 


sina 


this yields 

d 2 , , 


(cosa) 2 -^(j)i{z) + (1 + sina 1 ) 2 ^ 


(z = 


(31) 

-2TTpi(z) 

(32) 


and, consequently, 
d 2 

Qp [Mz) - <h(z )] = —27r(cos a) 2 (pi(z) - p 2 (z)) . 

(33) 


This means, if an area of zeros of the partition function 
overlaps the real time axis, the second derivative of the 
free energy is discontinuous. 


on the real time axis, which are covered by areas of zeros. 
If the spectrum is gapped, the beginnings and end 
points of two consecutive intervals T* and T* +1 are 
equidistant with t b n+1 — t h n = and +1 — t e n = , 

where qt,/ e are momenta minimizing/maximizing E £ on 
the domain given by \B^\ = 1. The length of the single 
intervals T* increases linearly with n. However, if the 
spectrum is gapless, all those intervals extend to infinity. 

The condition (1351) allows for a physical interpretation; 
namely, the occurrence of DQPTs is through a continuity 
argument related to non-thermal mode occupation. - In 
BCS-type models the mode occupation is given by 

(ng) = (atdjj) = sin 2 (arctan.Bg) . (37) 


This means for all modes q , where the condition (1551) is 
satisfied, the mode occupation is (n$) = 1/2. Let us 
assume that for any two points in K there exists a path 
connecting both points, along which (nr) is continuous^ 
and the existence of modes with (ng) < 1/2. Both as¬ 
sumptions should be true for physically relevant mod¬ 
els and were found to hold in all cases considered in 
the Kitaev model. Then, we can set up the follow¬ 
ing chain of consequences: through the continuity con¬ 
dition, the existence of non-thermally occupied modes 
k + with (ng + ) > 1/2 implies the existence of modes q 
with ( ng) = 1/2. This in turn is equivalent to the fulfill¬ 
ing of the condition (1351) . which implies the occurrence 
of DQPTs in the time evolution. The mode occupation 
(ng + ) ^ 1/2 is non-thermal in the sense that it cannot 
be realized by Fermi-Dirac statistics with positive tem¬ 
perature. 

An equivalent formulation of condition (1351) is 


B. 2d BCS-type models 


In the Kitaev honeycomb model (and general BCS- 
type models) the partition function is given by (BUI) , i.e., 
the zeros in the complex plane are given by 


(k) = In + in(2n + 1) 


n e Z . (34) 


At this point it becomes obvious, that in the thermody¬ 
namic limit the double product over k x and k v in (l20l) 
leads to dense areas of zeros in the complex plane, since, 
generally, dk x z n (k) / ±dk y z n (k). These areas of zeros 

cover parts of the real time axis (z = it) if Re(z n (k)) = 0, 
i.e., if 


3 q G K such that = 1 . (35) 

Dubbing the B~ = 1 isoline £>i C K, there are intervals 


,* = (2n + 1 )tt 
n 2 E Bl 


(36) 


3g G K such that Aq-(a 0 )A^(ai) + e$(a 0 )e$(ai) = 0 , 

(38) 


where a is the quench parameter of the BCS-type Hamil¬ 
tonian, (see appendixlBl) . From this it becomes clear, that 
after quenching to a gapless phase there are zeros of the 
dynamical partition function on the real time axis, since 
E$(a) = 0 O £q(a) = A$(a) = 0. In the mode occu¬ 
pation picture this can be interpreted as follows: when 
quenching to a gapless phase, excitations cost no energy; 
thus, any quench produces inverted mode occupation. 

As discussed in the previous section, areas of zeros cov¬ 
ering the real time axis result in jumps in the second time 
derivative of the dynamical free energy density if there is 
a jump in the density of zeros. Eq. (1511) gives a “layer” 
of zeros for every n G Z. Therefore, our total density 
of zeros is a sum of the densities of the individual “lay¬ 
ers”, pz(z) = J2 n Pz( z )- Th e single layer densities are 
given the Jacobi determinant of the change of variables 
k —> z n (k)r— 


n G Z 
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P n Az) = 




dk a 
dlm(z n ) 
dk x 


dRe(z n ) 

dk y 

aim (z n ) 
dk v 


-1 


' ^ B \ 

2 E r Bl 
„ fc k 


HBZ) 

~2^t dk * E Z 

k 


* E l B l 

(2 n + l)7r 3 


d kErd k Bi - d kl 


(2 n + 1)tt 


2 El 

k 


^, E k I - 


EAB 


' d k Bl 

h y k 

2 E r Bl 
„ fc k 


In (Bi) 

k ~dkE, 


AE'iBl 

_ k k 

(2 n + l)7r 3 


dk x 

dBl 

dk x 


dk v 

8b\ 

_ k_ 

dk v 


2 El 

k 

-1 


(2 n + l)7r 


2 El 

k 


OkE, 


k = k(z n ) 


-1 


(39) 


At this point a more technical view of the zeros of the 
partition function is useful: the zeros z n (k) correspond 
to intersections of the isolines 


((2n + l)7rRe(z n (fc))\ 
exp -=- , 

\ Im( 2 : n (/c)) J 

(2 n + l)7r 

2Im(z n (A:)) 


(40) 


in the momentum plane. This means that the den¬ 
sity of zeros p™{z) diverges at the boundary, since there 
V-Eg || V£?2. Thus, when approaching the boundary of 
an interval T* from the inside of the interval, the second 
time derivative of Re [/(£)] will diverge. 


IV. DYNAMICAL PHASE TRANSITIONS IN 
THE KITAEV HONEYCOMB MODEL 

A. Zeros of the dynamical partition function in the 
Kitaev model 

In the Kitaev model not only quenches to the massless 
phase create inverted mode occupation. Also quenches 
across phase boundaries with final parameter J\ in a 
massive phase induce critical points in the real time 
evolution. It is physically reasonable to assume, that 
the mode occupation number (ng) is sufficiently well be¬ 
haved, namely, that for any two k 0 , k\ there exists a path 
7 : [0,1] -A [—7T, 7r] 2 with y(0) = ko and y(l) = k\ such 
that (ny( s )), s G [0,1], is continuous. We found this to 
be true for all considered cases. Under this prerequisite, 
the existence of a fully occupied mode fc + , ( n £ + ) = 1, 
implies that (rig) = 1/2 somewhere, because ( n k=0 ) = 0. 
(ri£ + ) = 1 corresponds to \Bj: + \ = oo and this happens 
when 

0 = uj: + (J 0 }u^ + (J 1 ) +vj: + (J 0 )v^ + {J 1 ) 

A 0^u s+ (J o )^fc + (^i)-M £+ (Ji)u E+ (Jo) . 

One possibility to fulfill this is 

±1 = e £+W = £ fc+(^i) 

E fe - + (J 0 ) E£ + (J!) ‘ 


(41) 

(42) 


Now, consider a quench ending in the a;-phase (Jf > 
Jf+Jf) and k+ = (t r, 0). Then e £+ (Ji) = 2(Jf-Jf+Jf) 
and e^ + (Ji)/E^ + (Ji) = —1. We find that at this point 
both quenches, starting from another massive phase, 

JZ<JE + J5 =► -| ± 7^: = 1 > (43) 

■^k+ 

and from the massless phase, 

JS<Jg + JS => | ± 77T = 1 > (44) 

^k+ (Jo) 


lead to non-analytic behavior because (l42l) is fulfilled in 
both cases. The same can be shown for quenches ending 
in the other massive phases, only k + needs to be cho¬ 
sen appropriately. This shows that in the Kitaev model 
occupation inversion is produced by quenches within the 
massless phase or quenches crossing phase boundaries. 

Figure [4] displays locations of the zeros of the 
Losclimidt echo in the complex plane given by eq. (13411 
for two quenches, one within the A x phase and one from 
the A x phase to the massless phase. Both panels include 
a phase diagram with an arrow indicating the quench pa¬ 
rameters Jo —> J\. The numerical values for the param¬ 
eters for this figure and all following figures are listed in 
Tab. U in the appendix. The zeros do indeed form areas, 
which are restricted to the left half-plane for the quench 
within the massive phase but cover parts of the real time 
(imaginary z) axis when J 0 and J\ lie in different phases. 


B. Real time evolution 

On the real time axis the rate function of the 
Loschmidt echo C[t) = \Z{it)\ 2 reads 


r(t) = - 


2ir 2 


dk x dk y In 


yj 1 + 2 ES cos(2 Ej;t) 


Bt 

k 


'0 —TV 


1 + B l 

k 


( 45 ) 
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Im. time Re(z) Im. time Re( 2 ) 


FIG. 4. Distribution of zeros of the Loschmidt echo in the 
complex time plane computed according to eq. © for two 
different quenches. The zeros form areas in the complex plane, 
(a) Quench within one phase. The zeros are restricted to 
the left half plane and no DQPTs occur, (b) Quench to the 
massless phase. The zero areas overlap the real time axis (i.e., 
imaginary a axis) and DQPTs occur at the intersections of the 
boundaries of the single areas z n {k ) with the real time axis. 
Time is measured in units of ^2 a Jf. 


and the time derivative is 


k (46) 

Figure [5] shows the time evolution of the rate function 
and its time derivative for various quenches obtained by 
numerical evaluation of the corresponding integrals. The 
gray-shaded areas in the plots indicate the intervals T* 
(cf. eq. (13611 ) of vanishing partition function. The re¬ 
sults exhibit the properties expected from the previous 
considerations. The rate function is smooth for quenches 
within the gapped phases; however, nonanalyticities oc¬ 
cur when phase boundaries are crossed in a quench or 
after a quench within the gapless phase. The beginning 
and end points of the critical intervals are equidistant, re¬ 
spectively, and for quenches ending in the massless phase 
the intervals extend to infinity. As anticipated, nonana¬ 
lyticities only show up at the boundaries of the critical 
intervals. Moreover, the nonanalyticities emerge as dis¬ 
continuities of 7(f), i.e., kinks in r(t). 

Note that the two plots in panels a) and c) show the 
time evolution of the rate function after the two quenches 
for which Fig. []] shows the location of the zeros of the 
partition function. 


r(t) = 


dfc T dfc, 


BiE % sm(2E s t) 


1 + 2 B?. cos(2 Ej:t) 


C. The Id limit 

In the limit J a —> 0 for any a £ {cc, y, z } the 2d Kitaev 
model dH) degenerates and becomes a set of separate Id 



time a Ji) t 


time ) t 


FIG. 5. Real time evolution of the rate function of the 
Loschmidt echo (H51i and its time derivative (Hilll for various 
quenches. Both were obtained by numerical evaluation of the 
integrals. The gray-shaded areas indicate sections of the real 
time axis that are covered by areas of vanishing Loschmidt 
echo (cf. Fig. |4|). If zeros of the Loschmidt echo cover parts 
of the time axis, also the density of zeros on the time axis, 
is included. Kinks in the time derivative of the rate 
function are observed when quenching across a phase bound¬ 
ary or within the massless phase. 


spin chains. Let us consider the case J z = 0. The vanish¬ 
ing of one of the other two parameters will give the same 
result due to the threefold symmetry. For J z = 0 the 
condition for nonanalyticities (l35l) is fulfilled at q with 


cos (q x - q y ) 


jx jx I jy jy 

J i ' ■’() J i 

JX jy JX TV ■ 
^0 J 1 ' J 0 


Along this line also the spectrum is constant, 


(47) 


E,J = 



+ Jf + 2Jf J\ 


jx jx _i_ jV jV 

• y ( ) ■’ 1 ' J {) l 

jx jV jx jV ’ 

J 0 "T J 0 


(48) 


Thus, the critical intervals T* defined in eq. (1M1) become 
critical points 


* _ (2 n + l)7r _ 2n + 1 
” " 2 = 2 


(49) 


on the real time axis. Figure [HI shows the rate func¬ 
tions for two different quenches with J§ = Jf = 0. 
The quench in Fig. (SJi does not cross a phase bound¬ 
ary and therefore the rate function is analytic. However, 
in Fig. 03 J 0 = (}; §, 0) and Ji = (§, |, 0) lie in differ¬ 
ent phases, and according to eq. flU) there are critical 















































































the ground state |</>o) to be non-degenerate, 



FIG. 6. Time evolution of the rate function r(t) of the 
Loschmidt echo and its time-derivative for two quenches in 
the effectively one-dimensional Kitaev model with J z = 0. 
The quench in (a) does not cross the phase boundary and the 
rate function is analytic. In (b) the phase boundary is crossed 
and discontinuities of f(t) occur at equidistant instances in 
time t„. 


times t* = 2 " +1 </§7r at which the rate function becomes 

A y o 

singular. In particular, these singularities are discontinu¬ 
ities in the first time derivative of the rate function, not 
in the second as in the genuinely two-dimensional cases. 
This observation underlines the fact that the continuity 
of the first derivative is inherent to the higher dimension¬ 
ality of the non-degenerate Kitaev model. 


D. The long time limit 

In a recent worki^ it was stated that if the Loschmidt 
echo C(t) supports analytic continuation, then the long 
real time limit Hindoo C(t) and the long imaginary time 
limit linv^oo C{—ir) coincide. In the context of dynam¬ 
ical phase transitions this gives rise to the conjecture, 
that the occurrence of DQPTs is closely related to the 
long time behaviour of C(t). Dynamical quantum phase 
transitions occur if the zeros of the dynamical partition 
function cross the real time axis. This means, r(z), the 
rate function of the Loschmidt echo, is non-analytic in 
the r > 0 half plane; therefore, the long imaginary-time 
limit and the long real time limit do not necessarily have 
to coincide. 

On the imaginary time axis z = t, with the eigenbasis 
of the quenched Hamiltonian \cj > n ), corresponding ener¬ 
gies E n , and c n = (4>„\ipi), 


C(-it) = \(4>i\e Ht \ ipi) 


2 


Ei 


\ 2 e~ En 


n 


2 


J^CnCn'e E "' T (</>„| </>„') 

71 , 71 / 


(50) 


Now, shifting the energy such that Eq = 0 and assuming 


lim £{-ir) = |c„| 4 = \(4>o\4’i)\ i = -F 4 , (51) 


where T = |(</>o|t/>i)| is the fidelity. Thereby, the 
Loschmidt echo is connected to the fidelity in the large 
imaginary time limit. 

As the Loschmidt echo is not well defined in the ther¬ 
modynamic limit, one should rather formulate eq. m 
in terms of the rate function: 

lim r(—ir ) = — lim lim — ln£(— it) 

t—>oo r —yoo N—too iv 

= — lim — lnJ -4 (52) 


According to the previous considerations this yields 

lim r(t ) = — lim -En.F 4 , (53) 

if r{z) is analytic in the r > 0 half plane. This is 
a quite remarkable result: In the long time limit the 
Loschmidt echo approaches a value given solely by the 
overlap of the initial state with the ground state of the 
post-quench Hamiltonian, although the stationary state 
will surely never be that ground state. Quenching in¬ 
evitably produces an excited state. Moreover, some in¬ 
formation about the initial state is preserved for all times. 

For the Kitaev model the fidelity is 


J-=|(^o>| = 

1 


n'fc (i+ B k a -k a %) io> 




o\i>o) 


d 

= exp l - 


N f cl k f 2 \ 

T / 4^2 ln V 1 + B z) 


Thus, if above conjecture is valid, we should find 

£^ r W = i / d2fcln ( 1 + jB l) 


(54) 


(55) 


for quenches within the massive phases. Figure [7] shows 
the long time behavior of the rate function for a quench 
within the A x phase and for a quench crossing phase 
boundaries; indeed, the rate function converges to the 
value given by the fidelity after the quench within the 
massive phase. In the other case, however, the rate func¬ 
tion seems to converge, but the value it approaches differs 
from the one given by the fidelity. Various other cases 
were checked and the behavior was always consistent with 
above mentioned conjecture. 

One can explain the convergence of the rate function 
heuristically based on the specific form given in eq. 63. 
The expressions for the long time limit of the rate func¬ 
tion in eq. (15511 and the definition of the rate function in 
eq. (1451) only differ in the nominators in the argument 
of the logarithm, which are 1 and 1 + B^ie~ 2lE z t , respec¬ 
tively. In the long-time limit the factor e -2 ®^ 4 oscillates 
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r(t) 


ir(t) — -iln^ 4 — 



0 5 10 15 20 25 30 

time (]T) q Ji)t 


FIG. 7. Long-time behavior of the rate function (a) after a 
quench within the A x phase and (b) after a quench crossing 
phase boundaries. The rate function converges in both cases. 
For the quench within the massive phase it indeed approaches 
the value predicted by Heyl’s and Vojta’s conjecture (1521) . and 
in the other case the limit lies well off that value. 


extremely fast as a function of k. If Bi is slowly changing 
compared to these oscillations and also small such that 

In (l + S|e“ 2iB ^ ss B^e~ 2iE ^ , (56) 


density^ But are these excitations sufficient to induce 
dynamical quantum phase transitions? 

Assume the parameter J(t) of the Hamiltonian is not 
quenched immediately from J(t < 0) = Jo to J(t > 
0) = Ji, but continuously according to some protocol 
with J(t < 0) = Jo and J(t > t r ) = J\. In this case, the 
state of the system must at any time still be of the form 
given in eq. (fT6l) . because at any time H(J(t )) can only 
excite both modes with opposite momenta fc, —k in one 
fc-sector. For t > t r the dynamics are the same as after 
a quench, however, the B? will depend on the details of 
the ramping protocol. As discussed above B £ is directly 
related to the mode occupation number (n^). Thus, in 
order to determine, whether DQPTs occur after ramping 
from Jo to Ji instead of quenching, it is sufficient to 
compute the mode occupation at t = t r . 

In order to get the mode occupation (n^) we make use 
of the fact that the total time evolution is simply made up 
by the time evolution of independent two-level systems in 
the single fc-sectors and the corresponding Hamiltonians 
are 


H^(J(t)) 


a g (j(i))\ 
2 W k (J(t)) e^J(t)) ) 


(57) 


In these terms the initial state is the ground state of 
H k(Jo), and the time 

evolved state | ipu(t)) can be obtained by numerical in¬ 
tegration of the Sclirodinger equation. The mode occu¬ 
pation number after the ramping is then given by the 
overlap 


( n k) = IW£IV>fc(ir))| 2 , (58) 


then the contributions of neighboring points in the mo¬ 
mentum plane will cancel in the integral and therefore 
both integrals eq. m and eq. H55K become equal. 
However, if the integrand is singular, there are areas 
where \B^\ ss 1 and therefore the contributions of close- 
by points do not necessarily cancel. As a result the values 
of the integrals differ. 

In the absence of DQPTs Eq. (15^1) can also be derived 
rigorously for BCS-type models by considering the Tay¬ 
lor expansion of the logarithm in the integrand of Eq. 
(1451) and computing the time averages of the single con¬ 
tributions in the power series. 


E. Ramping 


where = E^t). 

Figure [5] shows the occupation numbers for a quench 
within the gapless phase and for linear ramping with 

f Jo ,t < 0 

J(t ) = < Jo + (Ji — J 0 )t/t r , 0 <t <t r (59) 
^ Ji , t ^ 

and ramping period t r = 50. As expected from the pre¬ 
vious considerations the quench produces regions of non- 
thermally occupied modes in the Brillouin-zone. Such 
areas are also present after the ramping. This means, 
that also in the time evolution for times t > t r there will 
be dynamical quantum phase transitions. 


It is known that dynamical quantum phase transitions 
also occur if the Hamiltonian parameter is continuously 
ramped across a critical point instead of quenching it 
However, it is not clear what happens after ramping the 
parameter within a gapless phase. In a gapless phase 
the adiabatic theorem does not apply and it is known 
that also slow ramping can produce a non-zero defect 


F. Quenching the magnetic field 

In the presence of a magnetic field the phase B (cf. 
Fig. [2]) becomes gapped and at J x = J y = J z = J 
there exists a diagonal form of the Hamiltonian^ with 
spectrum E^(J 1 k) as given in eq. (1131) . that maps to the 
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Occupation numbers (n(k)) 



1 


0.5 


0 


FIG. 8. (a) Mode occupation (ng) after quenching from Jo 
to Ji within the gapless phase and (b) mode occupation after 
linearly ramping from Jo to Ji with ramping time t r = 50. 


general two-band form (see appendix I~Y 1 for details) 

H (J, «) = X) (%(«) • a) % , (60) 

k 


where a = (cr x ,a y ,cr z ) T is the vector of Pauli matrices, 

7 I = and 

1 /Re(A s (J,«))\ 

fe fc( J ’ K ) = 9 I Im(Ag(J, k)) I • (61) 

V W«) / 


The magnetic field /i is contained in the parameter k ~ 
hxh j v hz . It introduces topological order in the B-phase, 
characterized by the Chern number 


= Si 


br ■ (dh br x 9/, &,-) 
dfc s dfc. y v x fc - 


— 7T ./ —7T 


| 6 |= 


= sign(At) . 


(62) 


It was demonstrated that in such systems any quench 
crossing the boundary between topologically distinct 
phases induces dynamical quantum phase transitions . 7 

Figure [9] shows the time evolution of the rate function 
after quenching the magnetic field within one phase and 
between two topologically distinct phases. As expected 
the signature of dynamical quantum phase transitions 
shows up after the quench across the phase boundary. 


V. CONCLUSIONS 

We demonstrated how the domains of zeros in 2d BCS- 
type models differ qualitatively compared to Id systems; 
namely, in the thermodynamic limit, the zeros coalesce 
to dense areas of zeros in the complex time plane rather 
than lines. The covering of intervals of the real time axis 
by such areas of zeros indicates the existence of critical 



L T 


\ 

\ / 

.. 

b) 






01234567 80 12345678 

time 'T) t lirnr: J“) t 


FIG. 9. Time evolution of the rate function in the Ki- 
taev model with additional magnetic field after quenching 
the magnetic field, (a) Quench within one phase, ko = 
0.5 —>■ /ti = 0.1. (b) Quench across the phase boundary, 

Ko = 0.5 —>■ Kl = —0.1. 


points in the time evolution. We showed how this leads to 
dynamical quantum phase transitions as discontinuities 
in the second time derivative of the dynamical free energy 
as opposed to discontinuities in the first derivative known 
from Id systems. 

It was found that in the Kitaev honeycomb model dy¬ 
namical quantum phase transitions occur after quenches 
across the phase boundaries and after quenches within 
the gapless phase. It is to our knowledge the first time 
that DQPTs are found in a two-dimensional model after 
quenching without crossing an equilibrium phase bound¬ 
ary. In accordance with the general considerations re¬ 
garding the dynamics of BCS-type systems, DQPTs in 
the Kitaev model show up at the boundaries of the inter¬ 
vals T* on the time axis which are included in a domain 
of zeros of the partition function. At these points the 
DQPTs appear as kinks in the first time derivative of 
the rate function of the Loschmidt echo. As was shown 
to hold for any BCS-type model, the curvature of the 
free energy density diverges when the boundary of such 
an interval is approached from the inward. 

In degenerate versions of the Kitaev model, which ef¬ 
fectively constitute one-dimensional spin chains, already 
the first derivative of the rate function becomes discontin¬ 
uous after quenching across a phase boundary as known 
from other Id models^ This underlines the fact that the 
continuity of the first derivative is inherent to the higher 
dimensionality of the non-degenerate Kitaev model. 

Moreover, we found for the Kitaev model that, in ac¬ 
cordance with a conjecture concerning that matter^ the 
long time stationary state of the rate function of the 
Loschmidt echo has a close connection to the occurrence 
of DQPTs: if no DQPTs occur, i.e., if the rate function is 
analytic in at least a half of the complex plane, the rate 
function approaches a value given by the fidelity. If the 
rate function is, however, non-analytic, it does in general 
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not converge to this value. The fact that the long time 
limit in absence of DQPTs is given by the fidelity deserves 
particular notice, since the fidelity is the overlap of the 
initial state with the ground state of the quenched Hamil¬ 
tonian, but the approached stationary state is surely an 
excited state. 
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The examination of the mode occupation numbers af¬ 
ter ramping the parameter of the Hamiltonian instead 
of quenching it implies that the subsequent time evolu¬ 
tion also exhibits dynamical quantum phase transitions. 
Moreover, it was demonstrated that in the presence of 
a magnetic field quenches between the topologically or¬ 
dered phases induce dynamical quantum phase transi¬ 
tions, which was previously proven to be a general fea¬ 
ture in the dynamics of topologically ordered two-band 
models in Ref. @- 

DQPTs were to date not observed in experiments. The 
Loschmidt echo C(t), which shows the non-analytic be¬ 
havior, is not directly connected to a quantum mechan¬ 
ical observable. As mentioned above, the work density 
was suggested as an measurable quantity, which could 
show the signature of DQPTs. However, the Loschmidt 
echo becomes exponentially small with increasing system 
size, whereas, strictly speaking, DQPTs occur only in the 
thermodynamic limit. Experimentally measuring the sig¬ 
nature of DQPTs in work densities will therefore be very 
challenging. Experimental consequences of DQPTs in 
other quantities than the work distribution function are 
currently being investigated. It was shown that DQPTs 
are connected to timescales of the order parameter dy¬ 
namics in symmetry broken systems^ in particular the 
occurrence of DQPTs is directly related to the sudden 
transition from monotonous to oscillatory decay of the 
order parameter after a quench.' 24 - Alternatively, the pre¬ 
viously mentioned generalized expectation values could 
serve as measurable quantities^ Moreover, it was shown 
in a recent work that DQPTs in Ising spin models ex¬ 
hibit scaling and universality and numerical results in¬ 
dicate that signatures of the DQPTs can be found in 
the dynamics of spin correlations as power-law scaling, 
which is solely determined by the universality class <2^ 
This seems to provide a very promising opportunity for 
measuring DQPTs, since these quantities are accessible 
with current experimental techniques. 

The finding of DQPTs as discontinuities in higher or¬ 
der derivatives in higher dimensional systems raises the 
question of a classification of dynamical quantum phase 
transitions. Canovi et al£ suggested a formalism for such 
a classification. They related discontinuities in general¬ 
ized expectation values and coexisting solutions to a first 
order transition. Here, we found a discontinuity in the 
second derivative of the dynamical free energy density. 
In future work it should be investigated how the findings 
of discontinuities in higher derivatives of the dynamical 
free energy density tie in with their definition. 


Appendix A: Bogoliubov transformation and 
post-quench eigenbasis 


Consider the general Hamiltonian 

H ( a ) = Yl (kfc( a ) ' % ( A1 ) 

k 

with a the vector of Pauli matrices 



(A2) 


and 7 £ containing the creation and annihilation operators 


= 



(A3) 


For a = J and b%(J) = (A%( J)/2, 0, eg( J)/2) this gives 
the Hamiltonian of the Kitaev model as in eq. (f5| . The 
unitary transformation 



W (a)7 E 


( M a ) v n ( a ) ^ 

\~ v k( a T M fcO)*/ 



that brings the Hamiltonian (IA1I) into diagonal form is 
the Bogoliubov transformation. Plugging eq. m into 
eq. m and demanding that all off-diagonal terms van¬ 
ish yields 


0 =d z % {u % v_ % - u Zk v % ) 

+ ( d l - id \) u -k u k + ( d ! + id l) v -k v k > ( A5 ) 

where all the dependencies on a have been dropped for 
the sake of brevity. Since W is unitary, uj: and v % have 
the general form 

ug = cos ffge*** , vg = sin . (A6) 

By choosing 

h = ~4>U = ~\ arctan ( A7 ) 

the equation above becomes real, 

0 =d| ^sin cos 9 £ — cos sin O^J 

+ \d x y\ (coscos- sin0£sin0_£) , (A8) 
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where d xy = d~ — id y - was introduced. Then for a - = 
k k k —k 

(a_r)t to hold, 6t must be an odd function of k. Thus, 


0 = sin(20g) + \d x J\ cos(2(9g) 

14 " I 

=► tan(2%) = (A9) 


and this yields 


l“fc| 2 =cos 2 %= i 

kl 2 = sin2 % = \ 



(A10) 


In the end it is left to choose the signs appropriately such 
that ug = u_g and ug = -v_g, e.g., 


with 

t/g(a 0) cci) = Mg(ai)Mg(a 0 )* + fg( Ji)v^{a 0 )* 
Vg(a 0 ,ai) = Mg(a 0 )Mg(ai) - wg(ai)ug(a 0 ) ■ (A14) 

Since |7i) is the ground state of H(J 0 ) 

a g°IV4 = (^(ao.aijag 1 - Vg(a 0 , ai)a“ 1 - t ) |V>i) 

= 0 (A15) 

must hold. Moreover, the ground state of a BCS-type 
Hamiltonian has vanishing total momentum; thus, 

Hi) = ^ { l +- B fc( a o,ai)«g lta “ 1 fc t ) l°;«i) 

^Bg(a 0 ,ai)a2 lt a“ 1 - t j |0;ai) , 

(A16) 




» E = sig(All) 

With this transformation the Hamiltonian (1A1I) becomes 
diagonal 


ff (°) = E ^ (“s'“E - “W) < A12 > 

k 

with E % {a) = |6g(a)|. 

To compute the quench dynamics one needs the con¬ 
nection of the degrees of freedom a2° that diagonalise the 
initial Hamiltonian H (op) and the degrees of freedom a^ 1 
diagonalising the final Hamiltonian H(a.i). This connec¬ 
tion is given by two subsequent Bogoliubov transforma¬ 
tions, 



, f aZ° 

W{ ai )W{a o ) 1 [ 4t 

( U^{a 0 ,Q!i) Vg(ao,ai )\ f a T\ 
V-Vg(a 0 ,Q!i)* C/g(ao,ai)*y \74b7 

(A13) 


where the coefficients B g(cxo, an) are to be determined, N 
is a normalization constant, and |0; oq) denotes the vac¬ 
uum of the post-quench fermions: ag 1 10; oq) = 0. Plug¬ 
ging eq. (1A16D into eq. (IA15I) yields 

(u % (a 0 ,a{)af - Vg(a 0 ,ai)^) 

x II ( 1 + jB fe'( Q; o,ai)ag / lt a“ 1 g ) t ) |0) ttl 

k' 

= (C/g(a 0 ,Q!i)Bg(Q!o,ai) - Vg(a 0 ,ai)) a“V 

x n , ( i+ - B£,(a ° ,a i) a i' ita -fe' t ) |0;ai) 

k'^k 

= 0 , (A17) 

which holds for 


Bg(a 0 ) an) 


Vg(ao,ai) 

t^g(a 0 , ai) 

M fc( a o)^fc( a i) - M fc(ai)r , fc(ao) 
Mg(a 0 )Mg(ai) + fg(ao)ug(aq) 


(A18) 


Appendix B: Condition for real time zeros of the 
partition function 

The condition B% = 1 (cf. eq. (l35l) f can be rearranged 
as follows: We have 
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1 = 1^1 





c fc(jo) 1 , e fc(^) 
^Po)\| ^(Ji) 


+ sign(A £ (J 0 )A g (Ji)) 


c fc(^) x e fc(^i) 

W)\ £^i) 





1 + 


e fe( J i) 


E S (Jo)\ E % {J\) 




e fcW V , f eg(Ji) \ 

£*KW ww 


Plugging this into the second last line we find the addi- which allows to write 
tional condition 

sign(Ag(J 0 )eg(Jo)) = -sign(Ag(Ji)eg(Ji)) , (B2) 

_I 


(Bl) 


e fcW V , ( ffcW V A x sigp(eg(ji)Ag(Ji)) 
Ek(Jo)J \Ez(Ji)J sign(eg(J 0 )Ag(J 0 )) 


(B3) 


Plugging in £g( J) 2 = eg(J) 2 + Ag(J) 2 yields 
1 = |-Bg(Jo, </i)| ^ Ag(Jo)Ag(Ji) + eg(Jo)eg(Ji) = 0 . 


(B4) 

So, for the emergence of a non-analyticity at a given time 
t = t* we get two (simplified) conditions: 


0 = Ag(J 0 )Ag(J 1 ) + eg(J 0 )eg(J 1 ) 


(B5) 

(B6) 


Appendix C: Quench/ramping parameters 
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Figure 

p 

J V o 

J§ 

J f 

jy 

■n 

KO 

Kl 

m 

0.8 

0.1 

0.1 

0.6 

0.2 

0.2 

0 

0 

Eh) 

0.8 

0.1 

0.1 

0.4 

0.3 

0.3 

0 

0 


Eh) 

0.8 

0.1 

0.1 

0.6 

0.2 

0.2 

0 

0 

Eh) 

0.8 

0.1 

0.1 

0.2 

0.1 

0.7 

0 

0 

Eh) 

0.8 

0.1 

0.1 

0.4 

0.3 

0.3 

0 

0 

Eh) 

0.4 

0.3 

0.3 

0.3 

0.4 

0.3 

0 

0 


Eh) 

0.1 

0.9 

0.0 

0.4 

0.6 

0.0 

0 

0 

Eh) 

0.25 

0.75 

0.0 

0.75 

0.25 

0.0 

0 

0 

Eh) 

0.9 

0.05 

0.05 

0.6 

0.2 

0.2 

0 

0 

Eh) 

0.9 

0.05 

0.05 

0.1 

0.8 

0.1 

0 

0 

Eh) 

0.8 

0.1 

0.1 

0.6 

0.2 

0.2 

0 

0 

Eh) 

0.8 

0.1 

0.1 

0.4 

0.3 

0.3 

0 

0 

Eh) 

1/3 

1/3 

1/3 

1/3 

1/3 

1/3 

0.5 

0.1 

Eh) 

1/3 

1/3 

1/3 

1/3 

1/3 

1/3 

0.5 

-0.1 


TABLE I. Summary of the quench parameters used for the 
figures in the main text. 
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